# This function converts ASCII Grid files (or gzipped versions of them) in the 
# origPath to GeoTiffs (with LZW compression, internal tiling and overviews)
# in the destPath.
# origPath - directory containing ASCII grid files
# destPath - directory to store the resulting GeoTiff files
# pattern -  file path pattern to match to find the files to convert.
#			 "\\.gz$" for XXX.asc.gz
#			 "\\.asc$" for XXX.asc
#			 If the files end in .gz then they will automatically be uncompressed.
# origEPSG - Spatial reference EPSG code of the ASCII grid
# destEPSG - Spatial reference EPSG code of the final GeoTiff. If different to
#			 the origEPSG then the grid will be rewarped.
# noData - 	 No Data to forcefully set when performing a warp. This is to overcome
#			 an apparent bug in gdalwarp no respecting the nodata value of the input
#			 file. If no warp is perform this value is ignored.
asciiToGeotiff <- function(origPath, pattern="\\.gz$", destPath, 
	origEPSG="4326", destEPSG="4326",
	noData="-9999") {
	
	
	library("R.utils");	# gunzip
	library("tools"); 	# built in library. file_path_sans_ext

	# List all file ending in .gz
	files <- list.files(path=origPath, pattern = pattern)
	for (i in 1:length(files)) {
		
		# ==================================
		#  Test if output Tif already exists
		# ==================================
		origFilePath <- paste(origPath,"/",files[i], sep="")
		decompressionNeeded <- grep("\\.gz$", origFilePath) == 1
		if (decompressionNeeded) {
			# Reduce XXX.asc.gz to XXX then add .tif to end
			finalGeoTiffFile <- paste(destPath,"/",
				file_path_sans_ext(file_path_sans_ext(files[i])),".tif", sep="")
		} else {
			# Reduce XXX.asc to XXX the add .tif
			finalGeoTiffFile <- paste(destPath,"/",
				file_path_sans_ext(files[i]),".tif", sep="")
		}
		if (file.exists(finalGeoTiffFile)) {
			print(paste("Skipping ", files[i], " Output already exists: ", finalGeoTiffFile, sep=""))
			next
		}

		
		# Record a list of temporary files that should be deleted once the final
		# product has been created.
		tempFiles <- c()
		
		# ================
		#  Decompression
		# ================
		# Check if decompression is needed
		if (decompressionNeeded) {
			srcFile <- origFilePath
			# Assumes the files are named XXX.asc.gz
			destFile <- paste(destPath,"/",file_path_sans_ext(files[i]), sep="") 
			
			print(paste("gunzipping: from:",origFilePath, " To:", srcFile))
			
			gunzip(srcFile, destname=destFile, remove=FALSE)
			tempFiles <- append(tempFiles, destFile)
			# Set the name of the source file for the next stage of processing
			srcFile <- destFile
		} else {
			srcFile <- origFilePath
		}
		
		# ===================
		# GeoTiff Conversion
		# ===================
		destFile <- paste(file_path_sans_ext(srcFile),"_Temp.tif", sep="")
		# Convert the ASCII grid to a GeoTiff using GDAL
		command <- paste(
			'gdal_translate -of GTiff -a_srs EPSG:',origEPSG,' ',
			srcFile, " ", destFile, sep="")
		print(command)
		system(command)
		srcFile <- destFile

		
		# ================
		#  Warping
		# ================
		if (destEPSG != origEPSG) {
			# rename the converted GeoTiff so that the output of the Warp conversion corresponds
			# to the output product name.
			srcFileRenamed <- paste(file_path_sans_ext(srcFile),"_EPSG",origEPSG,".tif", sep="")
			file.rename(srcFile, srcFileRenamed)
			
			command <- paste(
				'gdalwarp -t_srs EPSG:',destEPSG,' -srcnodata -9999 -dstnodata -9999 -co "TILED=YES" -co "COMPRESS=LZW" ',
				srcFileRenamed, " ", destFile, sep="")
			print(command)
			system(command)
			
			# The pre warped tif is now redundant so mark it for removal.
			tempFiles <- append(tempFiles, srcFileRenamed)
			srcFile <- destFile
		} 
		
		command <- paste(
			'gdaladdo -r average ',
			srcFile, ' 2 4 8 16 32 64',sep="")
		print(command)
		system(command)
		
		# Remove temp files
		file.remove(tempFiles)
		
		file.rename(srcFile, finalGeoTiffFile)
		print(paste("Completed ", i, " of ",length(files), sep="")) 
		
	}
}